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I.    INTRODUCTION 

The  efforts  to  arrive  at  coherent  descriptions  of  the  interactions  between 
oceans  and  the  structures  inserted  therein  have  a  long  history.  The  past  two 
decades  have  seen  an  explosion  of  interest  in  the  broad  subject  of  ocean 
hydrodynamics.  As  a  result  of  this  activity,  there  is  arising  an  improved  and 
more  realistic  understanding  of  the  physical  characteristics  of  some  time- 
dependent  flows  about  bluff  bodies  and  their  mathematical  formulation.  On 
the  one  hand  attention  has  been  focussed  on  controlled  laboratory 
experiments  which  allow  for  the  understanding  of  the  separate  effects  of  the 
governing  and  influencing  parameters,  and  on  the  other  hand  on 
mathematical  and  numerical  methods  which  allow  for  the  nearly  exact 
solution  of  some  wave  loading  situations. 

The  hydrodynamic  loading  situations  which  are  well  understood  are 
those  which  do  not  involve  flow  separation.  Thus,  they  are  amenable  to 
nearly  exact  analytical  treatment.  These  concern  primarily  the  determination 
of  the  fluid  forces  on  large  objects  in  the  diffraction  regime  where  the 
characteristic  dimension  of  the  body  relative  to  the  wave  length  is  larger 
than  about  0.2.  The  use  of  various  numerical  techniques  is  sufficient  to 
predict  accurately  the  forces  and  moments  acting  on  the  body,  provided  that 
the  viscous  effects  and  the  effects  of  separation  for  bodies  with  sharp  edges 
are  ignored  as  secondary. 

The  understanding  of  the  fluid-structure  interactions  which  involve 
extensive  flow  separation  and  dependence  on  numerous  parameters  such  as 


Reynolds  number,  Keulegan-Carpenter  number,  relative  roughness,  relative 
motion  of  the  body,  proximity  effects,  hydroelastic  response,  etc.  is  far  from 
complete  (Sarpkaya  &  Isaacson  1981).  There  are  several  reasons  for  this. 
First,  although  the  physical  laws  governing  the  motion  (the  Navier-Stokes 
equations)  are  well  understood,  valid  approximations  necessary  for  numerical 
and  physical  model  studies  are  still  unknown.  Even  the  unidirectional 
steady  flow  about  a  bluff  body  remains  theoretically  unresolved.  Much  of 
our  understanding  of  vortex  shedding  behind  bluff  bodies  came  from  steady- 
flow  experiments,  highly  idealized  models,  and  limited  numerical  solutions. 
Most  of  the  numerical  studies  based  on  the  use  of  the  Navier-Stokes 
equations  and  some  suitable  spatial  and  temporal  differencing  schemes  are 
limited,  out  of  necessity,  to  low  Reynolds  number  flows.  A  second  reason 
why  progress  has  been  slow  is  that  the  bluff  body  problems  involving  wake 
return  are  an  order  of  magnitude  more  complex  and  there  has  been  only  a 
handful  of  limited  applications  of  the  methods  based  on  Navier-Stokes 
equations. 

The  formation  of  a  wake  gives  rise  not  only  to  a  form  drag,  as  it  would 
be  the  case  if  the  motion  were  steady,  but  also  to  significant  changes  in  the 
inertial  forces.  The  velocity-dependent  form  drag  is  not  the  same  as  that  for 
the  steady  flow  of  a  viscous  fluid,  and  the  acceleration-dependent  inertial 
resistance  is  not  the  same  as  that  for  an  unseparated  unsteady  flow  of  an 
inviscid  fluid.  In  other  words,  the  drag  and  inertial  forces  are 
interdependent  as  well  as  time-dependent.  These  effects  are  further 
compounded  by  the  diffusion  and  decay  of  vortices  and  by  the  three- 
dimensional  nature  of  vorticity  due  to  turbulent  mixing,  finite  spanwise 
coherence,  and  the  random  nature  of  the  vortices  (which  give  rise  to  cycle-to- 


cycle  variations  and  numerous  flow  modes  even  under  controlled  laboratory 
conditions).  The  stronger  and  better  correlated  the  returning  vortices,  the 
sharper  and  more  pronounced  the  changes  are  in  pressure  distribution  on 
the  body  and  in  the  integrated  quantities  such  as  the  lift,  drag,  and  inertia 
coefficients. 

The  numerical  prediction  of  the  fluid-structure  interaction,  through  the 
use  of  finite-difference,  finite-element,  and  discrete-vortex  methods,  has 
attracted  considerable  attention  during  the  past  two  decades  and  produced 
laminar  flows  difficult  to  measure  and  turbulent  flows  hard  to  verify.  The 
reasons  for  this  are  relatively  simple.  Numerical  solutions  based  on  the  full 
Navier-Stokes  equations  are  not  stable  at  high  Reynolds  numbers  and  the 
instability  is  non-linearly  related  to  the  characteristics  of  the  base  flow,  the 
input  parameters,  and  the  discretization  conditions.  Also,  the  real  flow  at 
the  computed  Reynolds  numbers  may  be  turbulent,  at  least  in  some  regions 
of  the  flow,  and  the  numerical  experiment  does  not  imitate  the  physical 
experiment.  Furthermore,  the  observed  physical  and  numerical  instabilities 
do  not  necessarily  correspond  to  each  other.  Assuming  that  the  calculations 
for  a  given  flow  are  carried  out  at  sufficiently  small  Reynolds  numbers, 
where  the  flow  is  known  to  remain  stable  and  laminar,  one  quickly  discovers 
that  it  is  practically  impossible  to  measure,  to  any  credible  degree  of 
accuracy,  most  or  all  of  the  predicted  quantities  (except  the  Strouhal  number 
and  the  photographs  of  the  flow  patterns).  One  may  also  raise  the  question 
as  to  whether  the  two-dimensional  numerical  calculations  could  or  should 
ever  be  compared  with  physical  experiments  attempting  to  mimic  two 
dimensionality  through  the  use  of  various  passive  or  active  devices  (e.g.,  end 
plates  on  cylinders). 


Evidently,  one's  view  of  the  state  of  the  numerical  modelling  depends  to 
a  large  extent  on  one's  objectives.  For  example,  if  the  objective  is  to  obtain 
some  approximate  answers  and  flow  kinematics,  one  might  be  perfectly 
satisfied  with  the  existing  codes.  If  the  objective  is  to  match  the  measured 
and  calculated  results  (e.g.,  lift  and  drag  coefficients),  one  might  achieve  the 
desired  objective  by  fine  tuning  a  number  of  model  parameters  (e.g.,  the 
order  of  approximation  of  the  velocity  and /or  vorticity  gradients,  particularly 
near  the  wall,  mesh  size,  time  step,  type  of  discretization,  outer  boundary, 
just  to  name  a  few).  If  one's  objectives  are  to  perform  numerical  experiments 
for  sake  of  numerical  experiments,  with  no  concern  with  the  compatibility  of 
the  numerical  and  experimental  results,  then  one  can  objectively  asses  the 
model  instead  of  attempting  to  attribute  to  it  artificial  powers  of  prediction. 

As  far  as  the  turbulent  flows  are  concerned,  some  or  all  of  the 
predictions  of  the  numerical  calculations  for  a  given  flow  depend  on  the 
closure  model  used.  Some  models  do  better  than  others  for  some  flows  and 
worse  than  others  for  other  flows.  No  model,  however  sophisticated,  has  a 
corner  on  the  numerical  market.  Evidently,  the  solution  of  unsteady 
turbulent  flows  at  sufficiently  high  Reynolds  numbers  will  have  far  reaching 
theoretical  and  practical  consequences.  As  noted  above,  this  is  not  yet  the 
case,  and  the  solutions  must  necessarily  be  confined  to  cases  where  the 
accurate  prediction  of  physical  experiments  is  not  the  real  objective. 
However,  it  is  hoped  that  even  the  approximate  solutions  will  have  enough 
information  to  elucidate  the  physics  of  the  phenomenon. 

It  is  clear  from  the  foregoing  that  the  objectives  of  the  present 
investigation  are  to  carry  out  extensive  numerical  experiments  through  the 
use  of  the  vorticity-stream  function  form  of  the  Navier-Stokes  equations  and 


their  finite  difference  form,  on  co-existing  flows  (sinusoidal  or  non-sinusoidal 
oscillation  plus  steady  mean  flow).  The  expectations  are  that  the  results  will 
point  out  the  strengths  and  weaknesses  of  the  code,  for  the  particular  type  of 
formulation  used,  explain  the  reasons  between  the  various  numerical 
predictions  of  the  same  problem,  and,  hopefully,  shed  some  light  on  the 
physics  of  flows  heretofore  uncalculated. 


II.  BACKGROUND  STUDIES 

A  finite  difference  analysis  of  the  Navier-Stokes  equations  for  a 
sinusoidally-oscillating  ambient  flow  about  a  circular  cylinder  at  K  (Keulegan- 
Carpenter  Number)  =  UmT/D  =  5  (Re  =  1000)  and  K  =  7  (Re  =  700)  has  been 
attempted  by  Baba  &  Miyata  (1987).  Their  results  have  shown  that  the 
calculations  can  be  carried  out  only  for  short  times  (less  than  two  cycles  of 
flow  oscillation)  with  a  non-super  computer.  Murashige,  Hinatsu  and 
Kinoshita  (1989)  have  used  a  similar  method  to  analyze  three  cases  (K  =  5,  7, 
and  10)  at  higher  Reynolds  numbers  around  104.  The  flow  was  perturbed 
by  artificial  means  to  trigger  an  asymmetry.  At  K  =  10,  a  transverse  vortex 
street  appeared,  in  agreement  with  experimental  observations.  The 
numerical  simulation  of  steady  flow  past  a  circular  cylinder  undergoing  in- 
line and /or  transverse  oscillations  through  the  use  of  two-dimensional 
unsteady  Navier-Stokes  equations  was  undertaken  by  Lecointe  et  al.  (1987) 
for  relatively  small  amplitudes  (A/D  =  0.13).  Justesen  (1991)  presented 
extensive  results  obtained  from  a  numerical  solution  of  a  vorticity-stream 
function  formulation  of  the  Navier-Stokes  equations  for  the  flow  around  a 
circular  cylinder  in  planar  oscillating  flow  at  small  Keulegan-Carpenter 
numbers  in  the  subcritical  Reynolds  number  range.  Justesen  introduced  a 
straining  parameter  "a"  in  order  to  better  resolve  the  large  gradients  near  the 
cylinder  surface.  This  is  in  addition  to  the  logarithmic  straining,  commonly 
used  as  part  of  the  transformations,  for  a  better  resolution  of  the  gradients 
near  the  body.     Evidently,  Justesen's  transformation  for  a  =  0  defaults  to  the 


logarithmic  straining.  However,  "a"  becomes  another  disposable  parameter, 
dependent  on  at  least  K  and  Re.  Justesen  had  to  choose  judiciously  the 
value  of  the  straining  parameter  for  each  K  in  order  achieve  drag  and  inertia 
coefficients  in  satisfactory  agreement  with  those  obtained  experimentally.  A 
systematic  numerical  variation  of  the  governing  parameters  for  an  arbitrary 
U(t)  is  extremely  difficult. 

The  in-line  oscillations  of  a  cylinder  in  uniform  flow  (or  the  sinusoidally 
oscillating  flow  with  a  steady  mean  flow)  has  been  the  subject  of  intense 
interest  in  recent  years  (see,  e.g.,  Sarpkaya  &  Isaacson,  1981  and  Sarpkaya  & 
Storm,  1985)  in  connection  with  the  understanding  of  the  behavior  of  hot- 
wire anemometers  and  the  fluid  loading  of  structures  subjected  to  currents, 
gusts  and  other  types  of  unsteady  flows.  The  biassing  of  the  shedding  of  the 
vortices  by  the  current  causes  profound  changes  in  both  the  drag  and  inertia 
coefficients,  relative  to  their  no-current  values.  The  mobile  separation  points 
undergo  large  excursions,  as  much  as  120  degrees  during  a  given  cycle  of 
oscillation  over  a  circular  cylinder  (Sarpkaya  and  Butterworth,  1992).  These 
effects  are  further  compounded  by  the  diffusion  and  decay  of  vortices  and 
by  the  three-dimensional  nature  of  vorticity  due  to  turbulent  mixing,  reduced 
spanwise  coherence,  mutual-induction  instability,  and  the  random  nature  of 
vortices  which  give  rise  to  cycle-to-cycle  variations  and  numerous  flow 
modes  even  under  controlled  laboratory  conditions.  It  is  because  of  these 
reasons  that  the  present  work  is  dedicated  to  the  numerical  analysis  of  such 
flows. 


III.     NUMERICAL  REPRESENTATION 

A.    COMPUTATIONAL  METHOD 

Here  only  a  brief  description  of  the  computational  method  is  presented. 
A  more  in  depth  description  is  given  by  Wang  (1989)  and  Putzig  (1991). 

The  fluid  is  assumed  to  be  two-dimensional,  incompressible  and  vis- 
cous. The  governing  equations  for  the  solution  are  the  Navier-Stokes  equa- 
tions, with  the  stream  function  and  the  vorticity  as  independent  variables. 
To  achieve  a  higher  density  of  mesh  points  near  the  cylinder  surface,  the 
computational  domain  is  transformed  from  the  physical  plane  (polar  coordi- 
nates, see  Figure  1)  to  a  rectangular  plane  (Figure  2).  In  the  rectangular 
plane,  the  mesh  is  maintained  at  a  uniform  grid  spacing.  It  is  necessary  to 
have  more  mesh  points  closer  to  the  cylinder  surface  because  in  this  region 
the  gradients  of  both  the  vorticity  and  the  stream  function  are  the  largest. 

A  third-order  in  time,  second-order  in  space,  three-level  predictor-correc- 
tor finite-difference  scheme  is  used  to  solve  the  vorticity-transport  equation. 
A  Fast  Poisson  Solver  based  on  the  High  Order  Difference  approximation 
with  Identity  Expansion  (HODIE)  and  the  Fast  Fourier  Transformation  (FFT) 
provided  by  the  IMSL  mathematics  library  is  used  to  solve  for  the  stream 
function. 

The  unsteady  Navier-Stokes  equations  in  the  polar  coordinates,  as  de- 
fined by  the  vorticity  transport  equation  and  the  vorticity /stream-function 
equation  are, 
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to  and  \p  are  the  vorticity  and  the  stream  function,  v  is  the  kinematic  viscosi- 
ty, t  is  the  time  and,  r  and  6  are  polar  coordinate  directions  (see  Figure  1). 
The  velocity  components  in  the  r  and  6  directions  are  defined  by 


u  =  - — -  and  v  =  — .  (4) 

rd6  dr  v  ' 


The  boundary  conditions  for  the  physical  problem  are: 
(1)  no  slip  and  zero  normal  velocity  on  the  surface  of  the  cylinder 


ip  =  —  on  r  =  R  (5) 

dr 


and  (2)  the  potential  flow  at  infinity  is  defined  as 


ip  =  U(r-_)sinB  (6) 

r 


and  lo  =  0  at  r  =  oc.    U  is  the  external  flow  and  R  is  the  radius  of  the  cylin- 
der. 

The  coordinate  transformations  required  to  go  from  the  physical  domain 
to  the  computational  domain  are: 

r  =  R  •  e (a*  }  and  8  =  an.  (7) 

where  R  is  the  radius  of  the  cylinder  and  'a'  is  a  transformation  parameter. 

The  transformation  of  the  non-dimensionalized  vorticity-stream  function 
equations  and  their  finite  difference  form  through  the  use  of  the  central  dif- 
ference approximation  for  vorticity  and  a  two-step,  three-level,  predictor-cor- 
rector scheme,  with  a  third  order  accuracy  in  time,  are  described  in  detail  in 
Wang  (1989),  in  Fredrickson  (1990)  and  in  Putzig  (1991)  and  will  not  be  re- 
peated here. 

B.     CALCULATION  OF  THE  FORCE  COEFFICIENTS 

The  in-line  and  transverse  force  coefficients  are  determined  from  the  com- 
bined contributions  of  the  shear  and  pressure  forces  acting  on  the  cylinder. 
The  viscous  forces  are  calculated  from  x    =  \iu>.    The  total  in-line  force  then 

reduces  to 

F,L  =  -f(ps  -cos(e)Rde)  -  J>a)sin(u))Rde  (8) 

0  0 

and  the  total  lift  force  as 
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FL  =  -f  (ps  •  sin(6)RdB  )  -  f  uxocos  (w)  Rd6  .  (9) 
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After  dividing  the  in-line  and  the  lift-force  equations  by    (0.5pU  D)   and  de- 


fining 


-       2(PS"PJ 

Ps   "    ^2—  M 


the  force  coefficients  reduce  to 
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CIL  =  ~\  f  (Pi  '  cos(6)Rd0  )  -_•  f  |uosin  (co)Rd6  (11) 
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and 
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CL  =  -\[(Vs  "sin  (e)Rde  )  -«-•  J>«cos(u))Rd6  .  (12) 

0  0 


The  pressure  coefficient  is  determined  from  the  Navier-Stokes  equations 
in  terms  of  dimensionless  vorticity.  Once  integrated  with  respect  to  6,  one 
has 
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Equation  (13)  is  substituted  into  equations  (11)  and  (12)  to  determine  the  nu- 
merical scheme  for  the  in-line  and  transverse  force  coefficients, 
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The  radial  derivative  of  the  vorticity  on  the  surface  of  the  cylinder,  ap- 
pearing in  Equations  (14)  and  (15),  is  determined  through  the  use  of  discrete 
pointwise  approximations  of  various  orders,  ranging  from  second  to  tenth  or- 
der.   For  a  second  order  approximation,  one  has 
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For  higher  order  polynomials  Equation  (16)  may  be  written  as 
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+  0(A£n)         (17) 


in  which  the  coefficients  A  -  K  are  given  in  Table  1  below. 
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TABLE  1:  COEFFICIENTS  OF  THE  POLYNOMIAL  IN  EQ.  (17) 


n=2 

n=4 

n=6 

n=8 

n=10 

A 

-3/2 

-25/12 

-49/20 

-761/280 

-7381/25! 

B 

2 

4 

6 

8 

10 

C 

-1/2 

-3 

-15/2 

-14 

-45/2 

D 

4/3 

20/3 

56/3 

40 

E 

-1/4 

-15/4 

-35/2 

-105/2 

F 

6/5 

56/5 

252/5 

G 

-1/6 

-14/3 

-35 

H 

8/7 

120/7 

I 

-1/8 

-45/8 

J 

10/9 

K 

-1/10 
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C     CALCULATION  OF  THE  DRAG  AND  INERTIA  COEFFICIENTS 

If  one  were  to  associate  the  total  force  with  a  velocity-square-dependent 
drag  force  and  an  acceleration-dependent  inertial  force  then  the  coefficient  as- 
sociated with  the  latter  may  be  interpreted  as  some  measure  of  the  added 
mass.  But  one  must  bare  in  mind  that  such  a  decomposition  is  far  from  be- 
ing unique. 

It  has  been  customary  to  express  the  fluid  force  acting  on  a  body  mov- 
ing in  a  fluid  otherwise  at  rest  as 


F(t)  =  iPC^Ap|{Uo  +  U(t)}|{Uo  +  U(t)}  +Pk|;vikj(t)  (18) 


where  U0  represents  the  steady  velocity;     U(t),     the     time-dependent     oscilla- 
tions;   C)j    ,  the  Fourier-averaged  drag  coefficient  and  ku.    ,  the  Fourier-aver- 

U  1  J 

aged  added-mass  coefficient.    It  is  customary  to  use  an  inertia  coefficient  Cm 

for  a  fluid  in  motion  about  a  body  at  rest  through  the  use  of  C     =  1  +  k". 

The  Fourier  averages  of  the  drag  and  added-mass  coefficients  over  a  peri- 
od of  T  may  be  calculated  by  multiplying  both  sides  of  Equation  (8)  once 
with  U(t)  and  once  with  dU/dt  to  yield 

T 

2JT(t)U(t)dt 

C^  -  T 2 (19) 

pApJ|{U0  +  U(t)}|{U0  +  U(t)}U(t)dt 


and 
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T 

JT(t)djU(t)dt 
(Cm  -  1)  -  kj,  -  -i, (20) 

pVf[iLuo)]2d, 

0 


which  may  be  evaluated  readily,  provided  that  sufficiently  reliable  data  are 
available  for  F(t),  U„,  U(t),  and  dU(t)/dt. 

A   simple   dimensional   analysis   of  the   flow   under  consideration  shows 

that  the  time-averaged  force  coefficients  (C"  and  ku.)  are  functions  of  a  rela- 
tive amplitude  or  Keulegan-Carpenter  number,  Mach  number,  Reynolds  num- 
ber,   and    a    parameter    involving    U0    (e.g.,    U    (T/D)    or    UQ/[U(t)]        ). 

There  are  numerous  possibilities  regarding  the  definitions  of  the  relative  am- 
plitude or  Keulegan-Carpenter  number  and  the  Reynolds  number.  The  pur- 
pose of  the  search  for  a  more  suitable  Keulegan-Carpen+er  number  and /or 
Reynolds  number  is  to  enhance  the  correlation  of  the  data  to  reduce  the  num- 
ber of  the  governing  parameters,  possibly  eliminating  U0(T/D)  as  an  inde- 
pendent parameter.  The  list  of  possible  Reynolds  numbers  and  Keulegan- 
Carpenter  numbers  is  long  and  will  not  be  given  here.  Suffice  it  to  note  that 
the  two  force-coefficients  for  the  flow  about  a  cylinder  may  be  written  as 

/Cu 
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or  as 
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^   =  f.(K+,Re+,VK)  (22) 


'i 


in  which 
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Re+  =  Re  '  1  +      °  >  =  Re  +  -1-  .  (23) 


The  purpose  of  the  present  calculations  was  not  to  provide  a  detailed 
comparison  between  the  measured  and  calculated  forces  but  rather  to  at- 
tempt to  establish  a  relationship  between  the  shedding  of  vortices  and  the  rel- 
ative magnitude  of  the  current.  The  particular  values  of  K  and  Vr  chosen  for 
the  calculations  (K  =  4  to  6,  Vr  =  0.0  -  1.2)  was  one  for  which  some  experi- 
mental data  were  available  at  comparable  p,  Re,  and  Vr  values. 
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IV.    DISCUSSION  OF  RESULTS 

The  numerical  experiments  were  carried  out  through  the  use  of  a  VAX- 
2000,  a  VAX-3520,  a  CRAY  Supercomputer  and  the  IMSL  Mathematical 
Library.  The  solution  procedure  and  technique  have  been  validated  for 
several  types  of  unsteady  flows,  i.e.,  impulsively-started,  suddenly-stopped, 
and  uniformly-decelerated  flows  before  applying  it  to  co -existing  flows 
(oscillatory  flow  or  sinusoidally-oscillating  flow  with  mean  velocity). 
Excellent  agreement  with  flow  visualization  and  experimentally  determined 
drag  and  lift  coefficients  has  been  obtained  for  both  symmetric  and 
asymmetric  wake  solutions.  It  is  this  validation  that  led  to  the  exploration  of 
the  characteristics  of  sinusoidally-  and  non-sinusoidally-oscillating  flows 
superimposed  on  a  mean  velocity. 

The  sinusoidal  oscillation  was  specified  by    U  =  UQ+ Umsin  (2jtt/T)    in 

which  U0  is  the  steady  mean  velocity  and  Um  is  the  amplitude  of  sinusoidal 
oscillations.    The  non-sinusoidal  oscillation  was  represented  by 

U  =  Uo  +  Um[320166cos(B)  -  76230cos  (30)  +  22869cos  (56) 

-5445cos(7B)  +847cos(90)  -63cos(110)  ]/(262144)  (23) 

in  which  6  =  art  -  ji/2.  This  profile,  without  the  mean  velocity  U(y  was  first 
used  by  Justesen  and  Spalart  (1990)  in  connection  with  the  turbulence 
modeling  of  oscillatory  boundary  layers.  They  have  referred  to  it  as  the 
"steep"  velocity  variation  since  it  yields  a  quasi-steady  interval  and  an 
interval  with  a  very  large  pressure  gradient  in  each  half-cycle  (see  Figure  3). 
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As  noted  by  Jensen  and  Spalart,  the  coefficients  in  equation  (23)  were 
determined  by  requiring  that  the  first  ten  derivatives  be  zero  at  t  =  0. 

The  flow  was  perturbed  by  changing  the  direction  of  the  ambient  flow 
sinusoidally  (with  an  amplitude  of  one-half  of  a  degree)  during  the  first  cycle 
of  the  oscillation.  The  amplitude  of  the  sine  wave  was  the  only  free 
parameter.  It  is  worth  noting  that  this  type  of  disturbance  gradually  returns 
the  perturbed  quantity  to  its  initial  state. 

Numerical  experiments  have  been  carried  out  in  the  range  of  K  =  4  to  6, 
P  =  200,  Re  =  800  to  1200,  At  =  0.002  and  for  various  values  of 
V    =  U  /U    ,  as  shown  in  Table  2. 

r  o         m 

TABLE  2:  THE  RANGE  OF  THE  GOVERNING  PARAMETERS 

7  8 


vr  / 

K  = 

4 

5 

6 

0.0 

X 

0.20 

X* 

0.40 

X 

X* 

0.50 

X 

X 

0.60 

X* 

X 

X* 

0.65 

X 

0.70 

X* 

X* 

0.80 

X 

X 

1.00 

X 

1.20 

X 

X 


*The  calculations  have  been  carried  out  for  both  the  sinusoidal 
and  non-sinusoidal  oscillations  at  the  indicated  K  and  Vr  values. 
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Evidently,  this  is  a  rather  limited  exploration  of  a  highly  complex 
problem  and  requires  much  more  numerical  and  experimental  work.  The 
purpose  of  the  present  calculations  was  not  to  provide  a  detailed  comparison 
between  the  measured  and  calculated  forces  but  rather  to  attempt  to 
establish  a  relationship  between  the  shedding  of  vortices  and  the  relative 
magnitude  of  the  current. 

Figures  4  through  9  show  at  time  t/T  =  9.0,  for  K  =  4  and  Vr  =  0.4,  the 

in-line  force,  the  transverse  force,  a  close-up  view  of  the  streaklines,  the 
complete  streaklines,  the  streamlines,  and  the  vorticity  field,  respectively; 
Figures  10-15  show  the  same  plots,  for  K  =  4  and  Vr  =  0.5;  Figures  16-21 
the  same  plots,  for  K  =  4  and  Vr  =  0.6;  Figures  22-27  show  similar  plots,  for 
K  =  4  and  Vr  =  0.65;  Figures  28-33  show  the  corresponding  plots,  for  K  =  4 
and  Vr  =  0.7;  and  Figures  34-39  show  the  same  plots,  for  K  =  4  and  Vr  =  0.8. 
Figures  40-45,  46-51,  52-57,  and  58-63  show  the  plots  cited  above,  at  time 
t/T  =  9.0,  for  K  =  6  and  Vr  =  0.6,  0.7,  1.0,  and  1.2,  respectively.  Figures  64- 
69  and  70-75  show  the  force  and  vortex  plots  for  the  steep-sine  oscillations 
for  K  =  4  and  Vr  =  0.6  and  0.7,  respectively.  Finally,  Figures  76-81  and  82-86 
show  similar  plots  for  K  =  6  and  Vr  =  0.6  and  0.7,  respectively. 

The  streaklines  show  that  for  relatively  small  values  of  Vr  ,  the  wakes 
are  on  both  the  upstream  and  downstream  sides  of  the  cylinder  and  exhibit 
highly  irregular  and  often  highly  stretched  vortices.  As  Vr  increases,  the 
wake  begins  to  occupy  more  and  more,  the  downstream  side  of  the  cylinder 
(see  Figure  13).  The  most  interesting  feature  of  the  flow,  however,  is  not  the 
position  but  the  shape  of  the  wake.  Figure  13  shows  that  the  width  of  the 
wake   increases   and   the  vortices   begin   to   arrange   themselves   along   three 
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rows.  Figures  18-20  for  Vr  =  0.6,  Figures  24-26  for  Vr  =  0.65,  and  Figures  30- 
32  for  Vr  =  0.7  show  that  the  wake  is  comprised  of  three  rows  of 
heterostrophic  vortices.  They  differ  only  in  detail  from  one  Vr  to  another  in 
the  narrow  range  of  Vr  from  about  0.6  to  0.7.  At  lower  Vr  values  (see 
Figures  6-8  for  Vr  =  0.4  and  Figures  12-14  for  Vr  =  0.5),  where  the  sinusoidal 
oscillation  is  relatively  important,  the  shedding  of  the  vortex  couples  become 
more  and  more  alternating.  One  pair  goes  to  one  side  of  the  cylinder,  the 
next  pair  goes  to  the  central  street,  and  the  third  pair  goes  to  the  other  side 
of  the  street.  Then  the  events  repeat  themselves.  At  higher  Vr  values, 
however,  (see  Figures  37  for  Vr  =  0.8)  the  central  vortex  pairs  become  weaker 
and  fairly  stretched  out.  The  vortex  pairs  rotate  and  orient  themselves  as  if 
they  were  going  to  be  part  of  an  ordinary  Karman  vortex  street. 

The  lift  coefficient  plots  show  (see  Figures  17,  23,  29,  and  35)  the 
increasing  asymmetry  of  the  wake  with  Vr  for  Vr  larger  than  about  0.7.  At 
lower  Vr  values  (see  Figures  16  and  22)  the  in-line  force  coefficient  remains 
essentially  constant  and  periodic.  At  larger  Vr  values  (see  Figure  34), 
however,  the  in-line  force  becomes  increasingly  asymmetrical. 

Table  3  and  Figure  87  show  a  comparison  of  the  calculated  and 
experimental  drag  and  inertia  coefficients  for  K  =  4  for  representative  values 
of  Vr  .  The  computed  values  were  based  on  an  eighth  order  polynomial 
representation  of  the  radial  vorticity  gradient.  As  expected,  the  inertia 
coefficients  agree  extremely  well.  As  far  as  the  drag  coefficients  are 
concerned,  the  agreement  is  not  as  good,  but  certainly  better  than  expected 
in  view  of  the  fact  that  p  was  200  in  the  calculations  and  about  600  in  the 
experiments.    Nevertheless,  the  trend  of  the  data  is  well  predicted. 
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TABLE  3.  COMPARISON  OF  MEASURED  AND  CALCULATED  FORCE 
COEFFICIENTS    FOR  K  =  4    AND  VARIOUS  VALUES  OF  Vr 

Vr:  0.4  0.6  0.8  1.0 

C^(exp)  1.9  1.8  1.7  1.5 

Cu  (cal)  1.92  1.86  1.83  1.60 

m  v         ' 

C^(exp)  0.73  0.85  0.87  0.91 

Cj(cal)  0.52  0.72  0.86  0.85 

Figure  88  is  a  representative  sample  of  the  flow  visualization  obtained 
with  K  =  4  and  Vr  =  0.7  at  t/T  =  20.  The  symmetric  growth  and  motion  of 
the  vortices  are  clearly  visible.  This  figure  should  be  compared  with  Figures 
30.  The  similarity  is  rather  striking  in  spite  of  the  fact  that  p  was  200  in  the 
calculations  and  about  650  in  the  experiments. 

As  noted  previously,  the  calculations  were  carried  out  for  larger  values 
of  K  in  order  to  determine  as  to  whether  the  unusual  vortex  wake  seen  in 
Figures  30-32  is  a  consequence  of  a  unique  combination  of  K  and  Vr  values 
or  whether  it  would  occur  at  any  K  value  as  long  as  Vr  is  chosen 
judiciously.  The  results  have  shown  that  (see  Figures  40-63)  the  three-row 
vortex  street  occurs  only  for  K  =  4  in  a  suitable  range  of  Vr  values  (about  0.6- 
0.7).  For  example,  a  comparison  of  Figures  31  and  49  shows  that  the  wakes 
for  K  =  4  and  K  =  6  are  significantly  different  for  identical  values  of  Vr . 


21 


The  calculations  through  the  use  of  the  "steep-sine"  oscillation  and  a 
steady  current  have  shown  (see  Figures  64-86)  that  the  three-row  vortex 
wake  does  not  occur  and  that  the  wake  is  comprised  of  a  series  of  relatively 
active  and  inactive  vortical  regions.  The  inactive  regions  result  from  the 
periods  of  rapid  acceleration  and  the  active  regions  from  the  time  intervals 
during  which  the  velocity  is  nearly  steady  (see  Figure  3).  More  significantly, 
however,  figures  such  as  Figure  82  show  that  the  repetitive  but  non- 
sinusoidal  nature  of  the  flow  causes  dramatic  changes  in  the  in-line  (see 
Figures  46  and  82)  as  well  as  the  transverse  force.  Parts  of  the  cycle  become 
"inertia  dominated"  and  parts  of  the  cycle  become  "drag  dominated."  Had 
the  flow  been  non-repetitive  and  non-sinusoidal  (e.g.,  Gaussian),  the  drag 
and  inertia  dominated  regions  even  at  high  K  values  would  have  randomly 
occurred  throughout  the  history  of  the  motion.  It  is  because  of  this  reason 
that  it  has  not  been  easy  to  evaluate  the  ocean  data  or  possible  to  compare  it 
with  the  laboratory  data,  obtained  under  idealized  conditions.  It  does  not 
appear  that  this  situation  will  change  in  the  near  future. 
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V.  CONCLUSIONS 


The  investigation  reported  here  warranted  the  following  conclusions: 

1.  Even  the  higher  order  finite  difference  formulations  of  the  govern- 
ing equations  based  on  the  vorticity/stream-function  formulation  of  the  Navi- 
er  Stokes  equations  can  be  solved  for  only  relatively  small  Reynolds 
numbers.    This  is  primarily  due  to  stability  and  computer  constraints. 

2.  The  numerical  experiments  with  pulsating  flows  (oscillation  plus 
steady  mean  flow)  for  K  =  4  and  K  =  6  yielded  total  force  coefficients  in 
good  agreement  with  those  obtained  experimentally. 

3.  For  K  =  4  and  relative  current  velocities  of  about  0.6-0.7,  the  vorti- 
ces shed  nearly  symmetrically  at  each  cycle  and  gave  rise  to  a  most  unusual 
three-row  vortex  street,  where  each  row  is  comprised  of  a  pair  of  het- 
erostrophic  vortices.  For  relative  current  velocities  larger  than  about  one,  the 
vortex  wake  returned  to  the  asymmetric  mode,  as  is  encountered  in  a  regular 
Karman  vortex  street. 

4.  For  K  larger  than  about  4,  the  three-row  vortex  street  did  not  occur 
either  for  sinusoidally-  or  non-sinusoidally-oscillating  co-existing  flow,  show- 
ing that  the  three-row  vortex  wake  is  indeed  quite  unique.  It  remains  to  be 
seen  as  to  whether  it  occurs  at  much  larger  Reynolds  numbers  (larger  (3)  for 
K  =  4  and  V    =  0.6-0.7. 
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5.  Extensive  flow  visualization  studies  yielded  vortex  patterns  in  close 
agreement  with  those  predicted  numerically  at  the  corresponding  relative  cur- 
rent velocities. 

6.  The  foregoing  numerical  experiments  could  not  have  been  possible 
had  it  not  been  due  to  the  availability  of  a  VAX-3520  and  a  CRAY  supercom- 
puter. It  is  also  realized  that  calculations  at  higher  Reynolds  numbers  and 
for  larger  numbers  of  cycles  of  flow  oscillation  will  require  extremely  large 
CPU  times  even  on  a  supercomputer. 
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APPENDIX 


Figure  1.  Grid  in  the  Physical  Domain 
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Figure  2.  Grid  in  the  Computational  Domain 


Figure  3.  "Steep  Sine"  Curve  Compared  With  Sine  Curve 
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Figure  4.  Inline  Force  Coefficient,  K  =  4,  Re  =  800,  Vr  =  0.4 


Figure  5.  Transverse  Force  Coefficient,  K  =  4,  Re  =  800,  Vr  =  0.4 
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Figure  10.  Inline  Force  Coefficient,  K  =  4,  Re  =  800,  Vr  =  0.5 
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Figure  11.  Transverse  Force  Coefficient,  K  =  4,  Re  =  800,  Vr  =  0.5 
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Figure  16.  Inline  Force  Coefficient,  K  =  4,  Re  =  800,  Vr  =  0.6 


Figure  17.  Transverse  Force  Coefficient,  K  =  4,  Re  =  800,  Vr  =  0.6. 


36 


o 

a' 

ii 

H 

d 
ii 


jTV^  2 


c 

a") 


a 

X 

c 
U 

GC 
I — 

u 

3> 


37 


© 

II 
H 

■«-> 

so 

o 

ii 


00 

II 
II 


173 

g 

2 

CD 
On' 


0) 

J- 


38 


© 

a 

ii 

H 

o 
n 


ii 

01 

II 

0> 

C 


I* 

d 
ri 

(Ll 


s, 


39 


o 

ON 

II 
H 

O 

II 


o 
o 

00 


PC 


VI 

3 

o 

c 
o 

u 


r-i 
CM 


40 


Figure  22.  Inline  Force  Coefficient,  K  =  4,  Re  =  800,  Vr  =  0.65 


Figure  23.  Transverse  Force  Coefficient,  K  =  4,  Re  =  800,  Vr  =  0.65 
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Figure  28.  Inline  Force  Coefficient,  K  =  4,  Re  =  800,  Vr  =  0.7 
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Figure  29.  Transverse  Force  Coefficient,  K  =  4,  Re  =  800,  Vr  =  0.7 
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Figure  34.  Inline  Force  Coefficient,  K  =  4,  Re  =  800,  Vr  =  0.8 
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Figure  35.  Transverse  Force  Coefficient,  K  =  4,  Re  =  800,  Vr  =  0.8 


51 


52 


••'  \ 


7 


tr; 


o 

c-' 

II 

H 


\J 


>NV 


—> 


V 


V-. 


r 


!    If 


<J 


u 


V 


7-" 


]>s 


CO 

d 
ii 

> 

o 


& 
II 

0) 

C 

to 
d 

s- 

m 
u 


i 


AM 


53 


54 


55 


Figure  40.  Inline  Force  Coefficient,  K  =  6,  Re  =  1200,  Vr  =  0.6 
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Figure  41.  Transverse  Force  Coefficient,  K  =  6,  Re  =  1200,  Vr  =  0.6 


56 


.;«•• 


S 

*•** 


t 


^v 


>.cvv 


•*■•►... 


»w 


^-...  ^,..    .M-r**-* 


«•»'  " 


i  \* 


.1 


o 

o-' 

ii 

H 

o 


o 

o 


II 
C 


XI 

- 

4-1 


(Li 
> 

■D 
V 
r 

U 

r4 

v 

I 


57 


o 

On 
II 

H 
d 


o 
o 


P* 
SO" 

II 

^. 
c« 
QJ 

c 

s 

TO 
<L) 

Jh 

j— > 

tri 

s- 


58 


59 


60 


9.B 


Figure  46.  Inline  Force  Coefficient,  K  =  6,  Re  =  1200,  Vr  =  0.7 
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Figure  47.  Tranverse  Force  Coefficient,  K  =  6,  Re  =  1200,  Vr  =  0.7 
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Figure  52.  Inline  Force  Coefficient,  K  =  6,  Re  =  1200,  Vr  =  1.0 
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Figure  53.  Transverse  Force  Coefficient,  K  =  6,  Re  =  1200,  Vr  =  1.0 
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Figure  58.  Inline  Force  Coefficient,  K  =  6,  Re  =  1200,  Vr  =  1.2 
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Figure  59.  Transverse  Force  Coefficient,  K  =  6,  Re  =  1200,  Vr  =  1.2 
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Figure  64.  Inline  Force  Coefficient,  K  =  4,  Re  =  800,  Vr  =  0.6 
("Steep  Sine"  Oscillation) 
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Figure  65.  Transverse  Force  Coefficient,  K  =  4,  Re  =  800,  Vr  =  0.6 
("Steep  Sine"  Oscillation) 
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Figure  70.  Inline  Force  Coefficient,  K  =  4,  Re  =  800,  Vr  =  0.7 
("Steep  Sine"  Oscillation) 
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Figure  71.  Transverse  Force  Coefficient,  K  =  4,  Re  =  800,  Vr  =  0.7 
("Steep  Sine"  Oscillation) 
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Figure  76.  Inline  Force  Coefficient,  K  =  6,  Re  =  1200,  Vr  =  0.6 
("Steep  Sine"  Oscillation) 
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Figure  77.  Transverse  Force  Coefficient,  K  =  6,  Re  =  1200,  Vr  =  0.6 
("Steep  Sine"  Oscillation) 
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Figure  82.  Inline  Force  Coefficient,  K  =  6,  Re  =  1200,  Vr  =  0.7 
("Steep  Sine"  Oscillation) 
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Figure  83.  Transverse  Force  Coefficient,  K  =  6,  Re  =  1200,  Vr  =  0.7 
("Steep  Sine"  Oscillation) 
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Figure  87.    Cd  and  Cm  versus  Vr ,  K  =  4,  b  =  200,  n  =  8. 

Cd  numerical  (solid  circles),  Cm  numerical  (solid  squares) 

Cd  experimental  (open  circles),  Cm  experimental  (open  squares) 
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